import "babel-polyfill";

const {E, random, exp, log, sqrt, sin, tan} = Math;

export var uniform = function(a, b) {
    return a + (b - a) * random();
};

export var bernoulli = function(p = 0.5) {
    return random() < p;
};

export var binomial = function(n, p = 0.5) {
    var result = 0;
    for (let i = 0; i < n; ++i) {
        if (bernoulli(p)) {
            ++result;
        }
    }
    return result;
};

export var negativeBinomial = function(n, p = 0.5) {
    var result = 0;
    while (n) {
        if (bernoulli(p)) {
            --n;
        } else {
            ++result;
        }
    }
    return result;
};

export var geometric = (p = 0.5) => negativeBinomial(1, p);

export var poisson = function(lambda) {
    var result = -1;
    var x = 0;
    while (x >= -lambda) {
        ++result;
        x += log(random());
    }
    return result;
};

export var exponential = function(lambda) {
    return -log(random()) / lambda;
};

export var gamma = function(alpha, beta) {
    if (alpha === 1) {
        return -beta * log(random());
    } else if (alpha > 1) {
        let ainv = sqrt(2 * alpha - 1);
        let bbb = alpha - log(4);
        let ccc = alpha + ainv;
        while (true) {
            let u1 = random();
            let u2 = random();
            let v = log(u1 / (1 - u1)) / ainv;
            let x = alpha * exp(v);
            let z = u1 * u1 * u2;
            let r = bbb + ccc * v - x;
            if (r + 1 + log(4.5) - 4.5 * z >= 0 || r >= log(z)) {
                return x * beta;
            }
        }
    } else {
        let x;
        while (true) {
            let b = 1 + alpha / E;
            let p = b * random();
            if (p <= 1) {
                x = p ** (1 / alpha);
            } else {
                x = -log((b - p) / alpha);
            }
            let u = random();
            if (p > 1) {
                if (u <= x ** (alpha - 1)) {
                    break;
                }
            } else if (u <= exp(-x)) {
                break;
            }
        }
        return x * beta;
    }
};

export var beta = function(alpha, beta) {
    var x = gamma(alpha, 1);
    var y = gamma(beta, 1);
    return x / (x + y);
};

export var arcsin = function() {
    return sin(random() * PI / 2) ** 2;
};

export var weibull = function(alpha, beta) {
    return alpha * Math.exp(-log(random()), 1 / beta);
};

export var extremeValue = function(mu = 0, sigma = 1, xi = 0) {
    var result;
    if (xi === 0) {
        result = -log(-log(random));
    } else {
        result = (Math.exp(-log(random()), -xi) - 1) / xi;
    }
    return mu + result * sigma;
};

export var normal = function(mu = 0, sigma = 1) {
    var result = sin(2 * PI * random()) * sqrt(-2 * log(random()));
    return mu + result * sigma;
};

export var lognormal = function(mu = 0, sigma = 1) {
    return exp(normal(mu, sigma));
};

export var chiSquare = function(n = 1) {
    var result = 0;
    for (let i = 0; i < n; ++i) {
        result += normal() ** 2;
    }
    return result;
};

export var cauchy = function(x0, gamma) {
    var result = tan((random() - 0.5) * PI);
    return x0 + result * gamma;
};

export var fisherF = function(m = 1, n = 1) {
    return (chiSquare(m) / m) / (chiSquare(n) / n);
};

export var studentT = function(n = 1) {
    return normal() / sqrt(chiSquare(n) / n);
};
